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Abstract 


A  modal  control  scheme  was  developed  for  linear  sys¬ 
tems  that  are  time  periodic .  The  scheme  was  then  applied 
to  a  satellite  in  orbit  about  the  earth-moon  Lagrange 
point,  L3.  A  theoretical  control  gain  prediction  was  also 
devised  and  then  verified  by  numerical  results.  The  feasi 
bility  of  the  control  scheme  was  then  verified  by  computer 
simulation  by  controlling  the  actual  non-linear  motions  of 
the  satellite  with  the  linearly  derived  controller.  The 
satellite  control  costs  were  found  to  be  30  percent  less 
than  the  control  costs  of  a  satellite  in  earth-synchronous 
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MODAL  CONTROL  OF  A  SATELLITE 


IN  ORBIT  ABOUT  L3 


I.  Introduction 


Background 

Since  the  advent  of  spaceflight  in  the  late  1950s  and 
early  1960s,  space  scientists  and  mission  planners  have 
sought  a  solution  to  the  problem  of  maintaining  the  orbits 
of  sensor-bearing  satellites.  Coupled  with  the  orbital 
problem  is  the  desire  to  use  the  minimum  number  of  satel¬ 
lites  that  will  satisfactorily  accomplish  the  mission  ob¬ 
jectives.  In  the  early  years  of  spaceflight,  missions  were 
altitude  limited  by  an  immature  booster  technology,  but  the 
technology  has  now  matured  to  the  point  where  virtually  any 
earth  orbit  is  conceivable.  A  further  complication  for 
space  scientists  is  the  fact  that  many  space  operations 
such  as  communications  and  strategic  reconnaisarica  satel¬ 
lites,  require  or  desire  total  coverage  of  the  earth's  sur¬ 
face  by  a  particular  constellation  of  satellites.  Current¬ 
ly,  such  space  operations  are  conducted  by  satellites 
positioned  in  earth-synchronous  orbits,  but  these  satel¬ 
lites  do  not  provide  total  global  coverage.  In  addition, 
synchronous  orbits  are  unstable  and  require  the  satellites 
to  carry  some  sort  of  orbit  maintaining  device. 
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The  lack  of  total  global  coverage  is  probably  the 
biggest  drawback  to  using  synchronous  orbits,  particularly 
for  communications  and  reconnaisance  applications.  Since 
synchronous  orbits  are  equatorial  and  the  orbits  are  at  an 
altitude  of  35,780  KM  (Ref.  6:111-4),  placing  three  satel¬ 
lites  in  synchronous  orbits  120°  apart  from  each  other 
would  provide  global  coverage  except  for  the  regions  of  the 
north  and  south  poles.  The  strategic  implications  of  no 
surveillance  or  communication  in  the  North  Pole  region  are 
obvious. 

In  an  era  in  which  the  two  "superpowers "  are  develop¬ 
ing  a  space  defense  capability,  it  is  also  important  to 
consider  a  possible  satellite  interception.  At  synchronous 
altitude,  the  transfer  from  the  earth's  surface  to  the  tar¬ 
geted  satellite  requires  roughly  six  to  eight  hours,  depen¬ 
ding  on  the  launch  point  of  the  interceptor,  the  location 
of  the  target,  and  most  importantly  the  boost  capability  of 
the  interceptor.  This  intercept  time  is  extremely  short 
•  (from  a  defensive  point  of  view)  when  one  considers  that 
before  evasive  action  can  be  taken,  the  defensive  forces 
must  determine  that  the  observed  launch  is  in  fact  an  in¬ 
terceptor,  the  interceptor's  target  must  be  identified, 
and  evasive  maneuvers  (provided  the  capability  exists)  must 
be  planned  and  executed.  Fcr  future  missions  of  strategic 
importance,  then,  earth-synchronous  orbits  are  simply  not 
acceptable . 
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Attempts  to  preserve  the  orbits  of  satellites  have 
thus  far  been  limited  to  improvements  in  control  systems. 
The  requirements  to  carry  orbit  maintaining  devices  and 
orbital  stability  devices  (RCS,  momentum  wheels,  control 
moment  gyros,  etc.)  can  combine  to  restrict  the  service 
life  of  a  satellite  in  the  following  manner.  Since  all 
launch  vehicles  that  are  currently  available  are  thrust 

% 

limited,  they  are  therefore  payload  limited.  This  payload 
limitation  forces  a  satellite  designer  to  often  limit  his 
operational  capability.  Certainly  this  problem  could  be 
negated  if  all  satellites  were  launched  by  a  Saturn  V  boos¬ 
ter,  but  such  a  plan  is  neither  practical  nor  economical. 
From  an  economic  standpoint,  a  cheaper  booster  permits  a 
more  expensive  payload.  Due  to  the  weight  limitation,  it 
is  desirable  to  allocate  the  vast  majority  of  available 
payload  weight  to  the  sensors  and  support  equipment  re¬ 
quired  to  accomplish  the  satellite's  primary  mission.  When 
control  systems  are  required,  the  weight  allocated  to  the 
sensors  is  decreased  by  the  weight  of  the  control  systems. 
When  the  control  system  fuel  supply  is  depleted,  the  satel¬ 
lite's  service  life  is  likewise  finished.  Adding  more  fuel 
to  extend  the  service  life  adds  more  weight,  further  de¬ 
creasing  available  sensor  weight.  Currently,  all  of  the 
above  factors  combine  to  force  the  satellite  designer  to 
make  design  tradeoffs  of  service  life  versus  mission  capa¬ 
bility. 


If  the  service  life  could  be  extended  by  decreasing 
the  use  of  a  control  system,  the  designer's  job  becomes 
much  easie::.  This  leads  one  to  search  for  higher  altitude 
orbits,  since  the  closer  an  orbit  is  to  the  earth,  the  fas¬ 
ter  it  decays.  Research  in  the  area  of  such  long-term 
stable  orbits  gained  increased  impetus  in  the  1960s  with 
the  initial  manned  and  unmanned  explorations  of  space  and 
consequent  discussions  of  eventual  space  colonization . 
Lagrange,  in  1772,  mathematically  demonstrated  the  exis¬ 
tence  of  equilibrium  solutions  of  the  equations  of  motion 
for  an  object  located  in  the  Sun- Jupiter  system.  His 
theory  was  finally  observationally  verified  with  the  dis¬ 
covery  of  the  Trojan  group  of  asteroids  in  1906.  Together, 
the  Trojan  discovery  and  Lagrange's  mathematics  demonstra¬ 
ted  that  there  are  points  in  space  in  the  vicinity  of  two 
attracting  bodies,  where  the  forces  on  a  much  smaller  mass 
balance.  Assuming  the  mass  of  the  third  body  is  suffi¬ 
ciently  small  that  it  has  no  effect  on  the  motion  of  the 
-  attracting  bodies,  the  balanced  forces  cause  the  third  body 
to  be  motionless  relative  to  the  two  attracting  bodies.  A 
great  deal  of  research  has  been  accomplished  on  this  "re¬ 
stricted  three-body  problem"  (Ref.  9:231-232). 

Since  the  Trojan  discovery,  Laqrange's  theory  has 
been  applied  to  the  earth-moon  system  and  the  solution  of 
the  mathematics  has  yielded  the  five  points  of  equilibrium 
shown  in  Figure  1.  These  five  "Lagrange"  points  (desig¬ 
nated  Ll  through  L5)  represent  theoretical  points  (in  a 
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Lagrange  Points  in  Earth-Moon  System 


coordinate  system  that  rotates  with  the  earth-moon  line)  in 
space  where  an  object  with  zero  velocity  (relative  to  the 
earth  and  moon)  may  be  placed  and  the  object  will  remain  at 
that  point  (which  rotates  with  the  frame) ,  thus  providing  a 
long-term  stable  orbit.  Motions  in  the  vicinity  of  the 
"collinear  points",  LI,  L2 ,  and  L3,  are  inherently  unsta¬ 
ble,  since  the  eigenvalues  of  the  system  equations  for  the 
collinear  points  have  at  least  one  eigenvalue  with  a  posi¬ 
tive  real  part  (Kef.  7:427).  The  same  analysis  shows  that 
the  "triangular  points",  L4  and  L5 ,  are  inherently  stable. 

The  above  theory,  though,  is  based  on  the  restrictive 
assumptions  of  circular,  planar  orbits  and  does  not  include 
the  very  sizeable  perturbations  of  a  fourth  body:  the  sun. 
When  the  "restricted  four-body  problem"  (earth,  moon,  sun, 
and  satellite)  is  analyzed,  the  triangular  points  become 
unstable  and  the  instabilities  of  the  collinear  points  be¬ 
come  worse.  One  concludes,  then,  that  the  long-term  stable 
orbits  at  the  Lagrangian  points  do  not  exist  when  the  more 
realistic  four-body  problem  is  considered.  The  points  do 
provide,  however,  a  starting  reference. 

Referring  again  to  Figure  1,  note  that  the  points  L3, 
1.4,  and  L5  form  an  approximately  equilateral  triangle  about 
the  earth.  If  it  were  possible,  positioning  satellites  at 
these  three  points  would  provide  the  desired  globa  L  cover¬ 
age  of  the  earth.  Since  the  Lagrange  points  lie  in  the 
earth-moon  plane  and  not  the  earth-equatorial  plane  of  syn¬ 
chronous  orbits,  one  satellite  positioned  at  each  of  the 
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points  L3,  L4 ,  and  L5  would  "see"  the  North  Pole  at  all 
times . 

Problem  and  Scope 

Since  the  Lagrange  points  are  attractive  from  a  de¬ 
signer's  point  of  view,  investigation  of  the  control  cost 
of  stabilizing  these  unstable  points  is  warranted.  The 
yardstick  for  cost  effectiveness  must  be  whether  or  not  the 
control  cost  is  comparable  to  the  control  cost  incurred  at 
synchronous  altitude.  Certainly  other  factors  merit  con¬ 
sideration  (i.e.,  additional  cost  to  boost  to  a  higher  al¬ 
titude,  advantage  of  total  global  coverage,  etc.),  but  this 
study  is  limited  to  the  issue  of  control  cost. 

Since  L4  and  L5  were  stable  prior  to  the  introduction 
of  the  perturbative  force  of  the  sun,  tney  require  less 
energy  to  stabilize  than  L3,  which  grew  even  more  instable 
with  the  sun's  force  included.  If  L3  can  be  efficiently 
stabilized,  then  L4  and  L5  should  likewise  stabilize  quite 
easily.  This  study,  therefore,  seeks  to  stabilize  a  satel¬ 
lite  in  the  vicinity  of  L3. 

Wiesel  used  the  techniques  of  Chapter  II  of  this  re¬ 
port  to  construct  a  periodic  orbit  about  L3 ,  which  is  shown 
in  Figure  2.  The  x  and  y  distances  of  Figure  2  are  rela¬ 
tive  to  the  origin  and  axes  shown  in  Figure  1.  The  orbit 
was  constructed  using  a  "very  restricted  tour-body"  model, 
in  which  all  bodies  lie  in  a  plane  and  the  orbits  of  the 
earth  and  moon  are  circular.  Smith  (Ref.  8)  reconstructed 
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the  Wicscl  orbit,  found  the  orbit  to  be  unstable  (as  expec¬ 
ted),  and  attempted  to  apply  a  linear,  constant-gain,  feed¬ 
back  controller  to  the  satellite  displaced  from  the  refer¬ 
ence  periodic  orbit.  After  attempting  various  gain  magni¬ 
tudes  and  various  feedback  quantities,  Smith  concluded  that 
the  control  law  was  inadequate,  based  on  the  control  costs 
he  observed. 

Smith  recommended  that  follow-on  studies  use  a  more 
realistic  model  of  the  four  Lody  system  and  a  different 
control  law.  This  study  inple/nents  both  of  these  sugges¬ 
tions.  The  system  model  was  modified  in  a  manner  described 
in  the  following  section  and  the  control  law  was  changed 
to  modal  control,  also  defined  in  a  later  section  of  this 
report.  With  the  above  modifications,  a  demonstration  of 
the  validity  of  orbits  about  L3,  L4 ,  and  L5  for  missions 
requiring  a  long  service  life  and  total  global  coverage  of 
the  earth  was  initiated. 
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II .  Problem  Analysis 


Modeling  Assumptions 

Smith's  research  was  based  on  a  model  of  the  three 
attracting  bodies  (earth,  moon,  and  sun)  that  placed  the 
moon  in  a  circular  orbit  about  the  earth  with  an  orbital 
period  of  27.321661  days.  The  earth's  orbit  about  the  sun 
was  also  assumed  to  be  circular  with  a  period  ol  363.256365 
days.  The  model  used  in  this  study  was  changed  from  the 
Smith  model  by  placing  the  moon  in  an  orbit  about  che  earth 
that  is  periodic  but  not  circular,  as  constructed  by  Wiesel 
(Ref.  10) .  The  true  motion  of  the  moon  about  the  earth  is 
not  actually  periodic,  but  the  Wiesel  lunar  orbit  more 
closely  describes  the  actual  motion  of  the  moon.  The 
earth  orbit  about  the  sun  was  still  modelled  as  circular 
with  the  same  orbital  period  used  by  Smith.  Other  inclina¬ 
tion  and  eccentricity  effects  were  initially  ignored  in 
this  study,  thereby  restricting  the  motions  of  the  four 
bodies  to  a  plane. 

Coordinate  System 

The  derivation  of  the  four-body  problem  Hamiltonian, 
which  is  described  in  the  following  section,  is  accom¬ 
plished  with  respect  to  the  coordinate  system  shown  in 
Figure  3.  Research  by  Heppenheimer  (Ref.  5)  and  others  has 
shown  that  the  periodicity  of  the  orbit  about  L3  is  not 
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Fioure  3.  Geometry  of  the  Four-Body  Problem 


readily  apparent  in  inertial  or  pseudo-inertial  coordi¬ 
nates;  therefore,  the  coordinate  system  chosen  for  this 
study  rotates  with  an  angular  velocity  equal  to  the  average 
angular  velocity  of  the  moon,  and  the  center  of  the  system 
is  located  at  the  earth-moon  center  of  mass.  This  frame 
allows  easy  visualization  of  the  periodic  nature  of  the 
orbit  about  L3. 

Derivation  of  the  Hamiltonian 

To  properly  derive  the  dynamics  for  the  restricted 
four-body  problem,  the  three  attracting  bodies  were  assumed 
to  have  their  actual  masses  while  the  satellite  shown  in 
Figure  3  (to  be  later  positioned  on  the  L3  orbit)  was  as¬ 
sumed  to  be  massless  (as  in  the  restrictee  three-body 
problem) .  This  assumption  assures  that  the  attracting 
bodies  dictate  the  motion  of  the  satellite,  but  the  satel¬ 
lite  has  no  effect  on  the  motion  of  the  attracting  bodies. 

As  in  the  derivation  of  any  Hamiltonian,  the  kinetic 
and  potential  energies  were  first  constructed.  The  iner¬ 
tial  point  in  this  derivation  was  assumed  to  be  the  earth- 
moon-sun  center  of  mass,  which  allows  the  inertial  position 
vector  of  the  satellite  to  be  written 

Fsat  =  *  -  SP-  P  (1> 

where  M  =  mi  +  m2  +  m3  (sum  of  the  masses  of  the  earth, 
moon,  and  sun,  respectively)  and  the  scale  factor  j-jp  lo¬ 
cates  the  inertial  point  of  the  system.  Taking  the 
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inertial  time  derivative  of  the  position  vector  in  the  ro¬ 
tating  frame  yields 


U.  /  v  — 

d  t  rSA'f  -  VSAT 


=  rR  + 


(i.sxR) 


,  m  3 . r—  ,  —  , 
(-jyj-1)  P  +  WX(- 


-)P  (2) 


where  -j-|prp  +  wx(— -p-)p  is  the  transport  velocity  of  the 
earth-moon  center  of  mass.  The  symbol  to  denotes  the  angu¬ 
lar  velocity  vector  of  the  rotating  frame,  which  is  actu¬ 


ally  the  moon's  mean  inertial  angular  velocity  vector,  and 
the  superscript  r  represents  differentiation  in  the  rota¬ 
ting  frame.  Kinetic  energy  is  defined  as 


VSAT* 


But  since  the  satellite  is  assumed  to  be  massless,  the  mass 
is  divided  from  both  sides  to  yield  the  specific  kinetic 
energy,  T  . 

Ts  =  +  <«xR>*  (wxR)  +~Jp-  rp  •"jja-  rp 

+  (wx-ppp)  *  (wx-jjp-p)  +  2  [ rR*  (ujxR)  ] 

+  2(rR*-^-  rf)  +  2[r£-  (ux-gip)  ] 

+  2  [  (ooxR)  •  -~jp-  rp]  +  2  [  (idxR)  •  (ux-^-p)  ] 

+  2[-g*  rt-  (wx-^p)J  ;  (4) 


The  potential  energy,  V,  of  the  satellite  is  due  to 
the  gravitational  attraction  of  the  other  three  bodies  and 
has  the  form 
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t  h 

where  r.  c._  is  the  magnitude  of  radius  vector  ‘.rou  the  i 

1  bn  X 

body  to  the  satellite  and  G  is  the  universal  gravitational 
constant.  As  in  the  case  of  the  kinetic  energy,  we  con¬ 
struct  a  specific  potential  energy  by  dividing  both  sides 
of  (5)  by  mSAT. 


V 

s 


— Gm  3 
|R  -  p| 


Gmi 


j  —2— — r 
'mx+m2 


+  R  | 


R- 


Gm2 


mi  -  , 

- 7 - r 

mi+m2  1 


(6) 


Since  potential  energy  is  not  dependent  on  rate  terms  (time 
derivatives  of  the  coordinates),  the  momenta,  p^ ,  of  the 
Hamiltonian  can  be  written 


3T 


where  q^  is  the  time  derivative  of  the  i  coordinate. 

Coordinates  chosen  were  the  x,  y.  and  z  components  of  the 

satellite's  position  in  the  coordinate  frame  of  Figure  3 

(R  ,  R  ,  and  R  ).  In  this  coordinate  frame,  the  p.  are 
x  y  z  i 

actually  inertial  velocity  components,  a  distinction  that 
became  very  important  for  the  derivation  of  the  controller. 
The  full  Hamiltonian,  then,  is  given  by 

H  =  Z  p.$.  -  L  =  l  p .$■  -  T  +  V  (8) 

i=l  1  1  i=l  1  1 

With  algebraic  simplification  of  the  above  equation  after 
all  terms  were  included,  the  Hamiltonian  reduces  to 


H 


h(Px2  +  Py2  +  P22)  +  u)(pxRy  -  PyRx) 


+  TT't  (rp  -  up  )p  +  (rp  +  wp  )p  ] 
M  x  y  Ky  x  *y 

_ Gm  2 _ 


x 

Gm3 

R-pl 


_ P.EU 

|  m  2 
mi  +m2 


r+R  I 


i  - 

nii  +ni2 


(9) 
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Since  the  motion  of  the  earth  about  the  sun  was  assumed  to 


be 

planar , 

no  z 

component  of  p  exists 

in 

the  Hamiltonian. 

Construction  of 

the  Periodic  Orbit 

Equations 

of  motion  are  formed 

in 

Hamiltonian  systems 

by 

forming 

= 

3H 

9Pk 

• 

• 

• 

V 

<N 

t— 1 

II 

n 

(10) 

K  = 

-8H 

k  “ ’  1  ,  2  |  •  ■  mg 

n 

(ID 

where  n  represents  the  number  of  generalized  coordinates 
describing  the  system  (Ref.  7:94).  To  generalize,  these 
equations  can  be  rewritten  as 

x(t)  =  f (x(t) ,t)  (12) 

where  x  is  a  state  vector,  and  for  this  study  is  composed 
of  the  coordinates  and  momenta  of  the  system.  The  vector 
f  contains  the  partial  derivatives  described  by  (10)  and 
(11),  and  is  a  function  of  both  the  states  and  time. 

The  goal  of  a  search  for  a  periodic  orbit,  then,  is 
to  find  a  solution  to  (12)  that  closes  on  itself  in  a  fi¬ 
nite  amount  of  time.  The  time  from  the  initial  starting 
position  until  the  orbit  closes  on  itself  is  defined  as  an 
orbital  period  and  is  denoted  by  x.  The  motion  from  t  to 
2t  in  a  periodic  system  should  be  a  repeat  of  the  motion 
from  0  to  t .  In  the  restricted  problem  of  three  bodies, 
the  Hamiltonian  itself  is  a  constant  of  the  motion,  which 
results  in  a  continuous  family  of  orbits  about  the  five 
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equilibrium  points.  This  is  not  the  ease,  however,  in  the 
restricted  four-body  problem,  since  the  Hamiltonian  becomes 
an  explicit  function  of  time.  If  periodic  orbits  about  the 
equilibrium  points  still  exist  in  the  four-body  problem, 
one  would  expect  the  orbits  to  have  a  one  to  one  resonance 
with  the  motion  of  the  moon.  Unpublished  research  by 
Wiesel  has  confirmed  that  the  one  to  one  resonance  condi¬ 
tion  does  apply  to  orbits  about  the  earth-moon  Lagrange 
points . 

Assuming  the  existence  of  a  periodic  solution  to 
(12),  xQ(t),  over  a  given  time  interval,  a  small  displace¬ 
ment,  6x,  may  be  added  at  the  initial  time,  t  .  The  dis¬ 
placed  solution  can  then  be  expanded  in  a  Taylor  series 
about  the  solution  xQ(t)  to  yield 


x  ( t )  =  xQ(t)  + 


3x (t)  | 

3x (tQ)  xQ(t) 


6x(to) 


+  higher  order  terms  (13) 


Truncating  (13)  at  first  order  and  defining 


<5x(t)  -  x(t)  -  xQ(t) 
one  obtains 


<5x(t)  =  <t  (t,tQ)  <5x(t0) 


a  variational  equation  for  the  system  where 


♦tt'V 


^*(t)  | 

3x(t0)  xQ(t) 


(14) 


(15) 


(16) 
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The  matrix,  4> ,  is  termed  the  state  transition  matrix  and  is 
a  function  of  only  the  initial  and  final  times  being  con¬ 
sidered.  The  nxn  4>  matrix  has  the  properties 

<I>(t.to)  =  <Mt,ti)<Mti  ,tQ)  (17) 

and 

<Mto,tJ  =  1  (18) 

where  I  is  the  nxn  identity  matrix. 

For  the  variational  equation,  (15),  to  be  of  value, 
a  differential  equation  for  the  propagation  of  iSx(t)  must 
be  derived.  Referring  to  the  equations  of  motion,  (12), 
and  the  equation  for  6x(t),  (14)  one  can  form  the  equation 

•  • 

xQ(t)  +  6i(t)  =  f (xQ(t)  +  6x(t),t)  (19) 

This  equation  can  also  be  expanded  in  a  Taylor  series  about 
xQ(t)  to  yield 

x  (t)  +  6x(t)  =  f(x  (t),t)  +  - j?-  |  _  x(t) 

dX(t)  XQ(t) 

+  higher  order  terms  (20) 
Substituting  Equation  (12)  into  (20)  produces 

6x(t)  =  A ( t ) 6x (t )  (21) 

to  first  order,  where 

A(t)  =  - —  |  _  (22) 

3x(t) 

Differentiating  (15)  with  respect  to  time  yields 
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(23) 


6x(t)  =  <1>  (t ,  tQ)  6x  (tQ) 

since  <5x(tc)  is  a  constant.  Equating  (21)  and  (23)  one 
obtains 

<Mt,t  )6x(t  )  =  A  ( t)  6x  (t)  (24) 

o  o 

But  6x(t)  is  defined  by  (15),  and  substituting  for  6x(t) 
produces 

'I1  (  L  ,  tQ)  iSx  ( tQ)  =  A(t)'l>(t,tQ)6x(to)  (25) 

*(t,to)  ®  A(t)4>(t,to)  (26) 

Numerically  integrating  equations  (12)  and  <26)  si¬ 
multaneously  with  chosen  initial  conditions  for  x(i:q)  and 
initial  conditions  for  <J>  defined  by  (18),  one  obtains  in¬ 
formation  on  the  solution  x  (t)  and  its  local  variation. 

o 

The  initial  conditions  chosen  for  x(t  )  are  actually 

o 

an  educated  guess  at  the  initial  conditions  for  the  actual 
orbit.  Integrating  with  these  initial  conditions  over  a 
’period,  t,  one  obtains  xq(t)  and  <J>(x,0).  Since  x(tQ)  was 
merely  a  guess,  iteration  will  be  required  to  find  the 
actual  orbit  and  accompanying  initial  conditions.  The 
error  between  the  initial  conditions  and  the  value  of  the 
states  after  integrating  over  t  is  expressed  as 

Ax  =  xq(t)  -  xQ(0)  (27) 

If  the  initial  conditions  a::»i  varied  by  6x(0),  the  resul¬ 
tant  6x(t)  will  not  necessarily  produce  the  sought  periodic 
orbit;  however,  if  we  let 
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(28) 


6x(0)  =  6x(i)  -  Ax 

the  periodic  orbit  will  be  produced.  Recalling  (15)  once 
ogain  and  substituting  it  into  (28)  yields 

Ax  =  4>  (x  ,  0 )  Sx  (0 )  -  6x  (0 )  (29) 

Ax  =  [<I>(t,0)  -  I]6x(0)  (30) 

The  above  vector  equation  can  be  viewed  as  a  linear  equa¬ 
tion  for  tSx(O)  in  terms  of  the  known  quantities  Ax,  <!■(•(, 0), 
and  1.  If  the  higher  order  terms  of  (13)  and  (20)  were  in¬ 
cluded  in  the  derivation,  only  one  iteration  step  would  be 
required  for  convergence;  however,  the  linearization  of  the 
two  equations  forces  an  iteration  until  the  solution  con¬ 
verges  to  a  chosen  tolerance. 

Implementation  of  the  Theory 

Using  the  theory  of  the  previous  section,  a  new  peri¬ 
odic  orbit  about  L3  was  constructed.  Equations  (12)  and 
(26)  were  numerically  integrated  simultaneously  over  one 
orbital  period  and  the  Ax  was  computed.  Since  the  system 
model  was  modified  only  slightly  from  the  Smith  orbit,  his 
initial  conditions  were  used  to  start  the  iteration  (Ref  8) . 
The  Ax  was  used  in  the  correction  equation  (30)  to  adjust 
the  initial  conditions  and  the  calculations  were  repeated 
until  convergence.  The  final  convergence  error  chosen  for 
the  orbit  was  1  X  10~10,  which  is  also  the  accuracy  of  the 
integration  routine  supplied  by  Wiesel.  Once  the  orbit  was 


19 


found,  a  Fourier  series  (Ref.  2:60)  was  fit  to  each  of  the 
states  of  the  orbit  (which  were  calculated  at  evenly  spaced 
points  around  the  orbit) ,  with  a  resulting  nineteen  coeffi¬ 
cient  Fourier  series  representation  of  the  orbit.  The  pur¬ 
pose  of  the  Fourier  series  representation  was  to  easily 
store  a  reference  orbit  for  control  analysis  (described  in 
Chapter  IV  of  this  report) .  The  new  orbit  is  shown  in 
Figure  4.  Again,  the  x  and  y  positions  shown  in  Figure  4 
arc  relative  to  the  carth-mcon  center  of  mass,  previously 
defined  as  the  origin  of  the  coordinate  system  (Figure  3) . 
When  the  new  orbit  is  compared  to  the  Wiesel  orbit,  one  no¬ 
tices  very  little  difference  between  the  two.  As  previously 
stated,  the  only  difference  in  the  construction  of  the  two 
orbits  are  the  two  models  of  the  moon.  The  addition  of  the 
periodic  versus  the  circular  orbit  for  the  moon  resulted  in 
a  slightly  smaller  orbit  in  both  x  and  y  directions.  The 
initial  conditions  for  the  Wiesel  orbit  and  the  new  orbit 
are  shown  in  Table  I. 


TABLE  I.  L3  Orbit  Initial  Conditions 


Wiesel  Orbit  New  Orbit 


x  position 
y  position 
z  position 
x  momentum 
y  momentum 
z  momentum 


(km) 

(km) 

(km) 

(km/sec) 

(km/sec) 

(km/sec) 


-383,519.08 

0.0 

0.0 

0.0 

-30 . 825652 

0.0 


-383,537.27 

0.0 

0.0 

0.0 

-30.825537 

0.0 
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Y  POSITION  (KM)  *102 

-41.00  -31.00  -21.00  -11.00  -i.00  9.00  19.00  29.00  39.00 

i  1  4*1  i  t  •  .  . 


Figure  4.  New  Orbit  About  L3 


III.  Floquet  Solution  for  the  Periodic  .System 

Floquet  Theory 

Unlike  the  restricted  three-body  problem  Hamiltonian, 
the  four-body  problem  Hamiltonian  is  an  explicit  function 
of  time.  Consequently,  the  Hamiltonian  itself  is  no  longer 
a  constant  of  the  motion.  Certainly  energy  is  conserved 
(T  +  V  =  E  =  constant) ,  but  no  other  constants  of  the  mo¬ 
tion  exist  (Ref  7:416).  The  equations  of  motion  derived 
from  the  four-body  Hamiltonian  are  also  found  to  be  time- 
varying  but  periodic,  as  one  would  expect  since  the  problem 
at  hand  is  an  orbital  problem.  Repeating  (21) ,  the  varia¬ 
tional  equations  of  motion  may  be  written  in  general  state 
vector  form  as 

• 

<5x  (t )  =  A (t )  6x  ( t )  (31) 

The  periodicity  of  the  system  is  contained  in  the  coeffi¬ 
cient  matrix  A(t)  such  that 

A (t )  =  A ( t+T  )  (32) 

where  t  is  once  again  defined  as  the  period  of  the  system. 
Differential  equations  such  as  (31),  which  have  variable 
yet  periodic  coefficients  can  theoretically  be  reduced  to 
the  case  of  constant  coefficient,  linear  differential 
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equations,  by  the  theory  of  Floquet  (Ref  4:60).  i'loquet's 
theory  states  that  the  general  solution  to  (31)  is  of  the 
form 


—  —  x  +• 

6x(t)  =  l  6  .E  .  (t) eA] 

j=l  3  3 


where  3j  and  X  ^  are  constants  and  the  E ^  are  periodic  func¬ 
tions  of  time  with  the  same  orbital  period,  x.  Actually, 


(33)  represents  a  special  case  solution  where  the  X^  are 
all  multiplicity  one.  The  n  X.  can  be  arranged  in  a  matrix 
of  Jordan  canonical  form  as  shown  below  (Ref  7:267). 


Previous  work  by  Smith  and  Wieoel  had  verified  the  simple 
multiplicity  of  the  X^;  therefore,  the  assumption  of  the 
form  of  (33)  is  warranted.  The  constant  X_.  is  called  a 
characteristic  exponent  or  Poincare  exponent  and  the  quan¬ 
tity  e*jfc  is  termed  a  characteristic  multiplier  of  A(t) . 
What  follows  is  a  method  for  computing  and  X^.. 

Assuming  one  desires  to  excite  only  one  mode  of  the 
system  (31),  let  0^=0  for  j  /  i  and  3^  =  1  which,  from 
(33) ,  yields 

6x(t)  =  Ei(t)cAit  (35) 
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From  ( 1 5 ) , 

6x(t)  =  4>  ( t ,  0 )  <$x  (0 )  (36) 

Equating  (35)  and  (36),  one  obtains 

Ei(t)eAit  =  <!>  (t,0)  6x(0)  (37) 

=  Ht,0)E.(0)  (38) 

Since  (t)  has  already  been  defined  as  a  periodic  vector, 
E.  Cl  )  =  iii  (0)  ,  such  that 

6x(t)  =  E^(x)eAiT  (39) 

=  IbCOJe^i1  =  <Mx,  0)^(0)  (40) 

Rearranging  (40)  yields  the  eigenvalue  problem 

{<D(t,0)  -  (eAiT)l}Ei(0)  =  0  (41) 

<t>  ( t  ,  0 )  is  called  the  monodromy  matrix  and  its  eigenvalues 
X  T 

are  e  .  Defining  the  monodromy  matrix  eigenvalues  as  a^, 
the  Poincare  exponents  can  be  found  by 

A.  =  -  ln(a.)  (42) 

IT  1 

The  eigenvectors  of  <t>(t,0),  therefore,  are  the  vectors 
Ei (0)  . 

To  compute  the  complete  solution,  E^(t)  must  be  com¬ 
puted  over  an  entire  period.  Since  Eb  (t)  is  periodic,  de¬ 
termining  E^(t)  over  one  period  determines  the  vector  for 
all  time.  Recalling  (35)  and  (31),  one  can  write 
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(43) 


|t(Ei(t)cXit)  =  A(t)  (13i(L)eAit) 

Now  define  A(t)  as  the  square  matrix  composed  of  the  n 
E^(t)  partitioned  into  columns 

A  (t)  =  [E .  ( t)  {  E2  (t)  |  .  .  .  !  E  (t)  ]  (44) 

i  n 

Now,  by  the  form  of  (43)  and  the  above  definition,  one  can 
show  that 

A (t)  cJt  +  A (t) JeJt  =  A (t) A (t) eJt  (45) 

A (t)  +  A ( t) J  =  A(t) A (t)  (46) 

The  above  equation  provides  an  ordinary  differential  equa¬ 
tion  for  A.  with  initial  corditions  E^(0),  previously  de¬ 
fined  as  the  eigenvectors  of  <J>(t,0). 

One  can  also  define  a  new  variable  £(t)  such  that 

£<t)  =  A (t) eJt  (47) 

and  from  (45) 

£  (t)  =  A(t)Ut)  (48) 

The  initial  condition,  £(0),  for  this  equation  is  the  E^(0) 
arranged  in  columns  or  A(0).  This  equation  has  the  same 
form  as  the  differential  equation  for  4>,  (26),  which  ver- 
fies  that  £  is,  in  fact,  4>,  to  within  a  coordinate  rotation. 

Implementation  of  the  Theory 

After  finding  the  Fourier  series  representation  of 
the  orbit,  the  stability  of  the  orbit  was  investigated. 
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Using  the  theory  of  the  previous  section,  the  system  eigen¬ 
vectors  and  Poincare  exponents  were  found  from  the  eigen¬ 
value  problem  of  (41)  combined  with  (42).  The  monodromy 
matrix  was  produced  by  intec  railing  (26)  when  the  periodic 
orbit  was  found. 

Smith's  analysis  of  the  Wiesel  orbit  (Ref  8)  demon¬ 
strated  that  the  Wiesel  orbit  is  unstable,  and  the  analysis 
of  the  new  orbit  indicates  the  same  instability.  As  pre¬ 
dicted,  one  of  the  system's  Poincare  exponents  is  a  posi¬ 
tive  real  root  (conjugate  root  is  negative  and  real) ,  while 
the  other  roots  are  pure  complex  conjugates.  Table  II 
below  demonstrates  the  slight  differences  between  the  Poin¬ 
care  exponents  of  the  two  orbits. 

Table  II.  Poincare  Exponents  of  the  Two  Orbits  (i.-imaginary ) 


Mode 


Poincare  Exponent 
Wiesel  Orbit  New  Orbit 


Planar 

Mode _ 

Planar 

Mode _ 

Out  of 
Plane  Mode 


0.0  -  1.03981 
0.0  +  1 . 0398i 
2.3932  +  O.Oi 
•2.3932  +  0  ■  Oi 
0.0  -  1 . 1248i 
0.0  +  1  .  1248i 


0.0  -  1 . 0392i 
0.0  +  1 . 0392i 
2.3921  +  O.Oi 
-2.3921  +  O.Oi 
0.0  -  1 . 1247i 
0.0  +  1 . 1247i 


In  both  orbits,  the  unstable  exponent  of  the  system  is  de¬ 
signated  A 3  and  is  the  root  the  controller  must  change. 


IV.  Control  Theor 


Review  o£  Previous  Work 

Previous  work  by  Wiesel  and  Smith  demonstrated  that 
there  are  three  modes  of  oscillation  for  the  periodic  orbic 
about  L3.  Two  of  the  modes  lie  in  the  orbit  plane  while 
the  third  mode  is  out  of  the  plane.  Wiesel  and  Smith  also 
discovered  that  the  out-of-plane  mode  and  one  of  the  planar 
modes  were  stable  yet  the  other  planar  mode  was  unstable, 
thereby  necessitating  some  form  of  control  on  the  satellite 
in  orbit  about  L3.  Smith  attempted  feedback  control,  at¬ 
tempting  both  rate  and  position  feedback  laws.  This  scheme 
was  unsuccessful  due  to  the  fact  that  all  modes  were  ex¬ 
cited  by  the  control,  which  is  contrary  to  the  desired  ef¬ 
fect  of  the  control  affecting  only  the  unstable  mode. 
Exciting  all  modes  resulted  in  control  costs  that  were  sig¬ 
nificantly  higher  than  the  costs  of  satellites  in  earth- 
synchronou;,  orbits.  Obviouuly,  ;hen,  an  alternate  control 
scheme  that  does  not  excite  the  two  stable  modes  must  be 
implemented . 

Modal  Control 

The  basic  premise  of  modal  control  is  to  perform  con¬ 
trol  on  a  preselected  number  of  system  modes  without  ex¬ 
citing  the  remaining  modes  of  the  system.  In  terms  of 
classical  control  theory,  modal  control  can  be  viewed  as  a 
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pole  placement  technique.  The  instabilities  of  the  system 
caused  by  roots  lying  in  the  right-half  of  the  complex 
plane,  are  removed  via  the  controller  moving  the  unstable 
roots  from  the  right-half  to  the  left-half  of  the  complex 
plane.  Simultaneously,  the  other  system  roots  remain  in 
their  original  locations  in  the  complex  plane. 

Since  the  Poincare  exponents  of  Hamiltonian  systems 
are  always  conjugate  pairs,  the  one  unstable  mode  of  the 
four-body  problem  implies  a  pair  of  roots  on  the  real 
axis:  one  root  in  the  right-half -plane ,  the  conjugate  root 

in  the  left-half-plane.  The  goal,  then,  of  the  modal  con¬ 
trol  applied  to  this  problem  was  to  move  the  one  root  out 
of  the  right-half-plane  while  the  other  five  roots  remained 
at  their  original  locations  in  the  complex  plane.  The  de¬ 
rivation  of  the  controller  follows. 

To  apply  control  to  the  system,  (31)  must  be  augmen¬ 
ted  to  yield  the  state  equation 

6x(t)  =  A  ( t)  <5x  ( t )  +  B  (t)u(t)  (49) 

where  u(t)  is  the  control  vector  and  B(t)  determines  to 
what  states  control  is  applied,  i.e.  B(t)  represents  the 
physical  reality  of  the  control.  To  insure  that  any  con¬ 
trol  applied  to  the  system  affects  only  the  mode  of  inter¬ 
est,  the  system  must  be  transformed  to  a  diagonal  set  of 
equations.  To  accomplish  this  objective,  define  a  new 
variable 
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n (t)  =  A"1 (t) 6x(t)  (50) 

where  r\  is  termed  the  modal  vector. 

Taking  the  time  derivative  of  the  above  equation  yields 

6x(t)  =  A(t)n(t)  +  A(t)n(t)  (51) 

Equating  (49)  and  the  above  equation 

A  ( t )  n"  ( t )  +  A(t)q(t)  =  A  (t)  6x  (t)  +  B(t)u(t)  (52) 

=  A(t)  A  (t)  n  (t)  +  B(t)u(t)  (53) 

Solving  for  n ,  since  a  differential  equation  for  propaga¬ 
tion  of  the  modes  is  desired,  yields 

n (t)  =  A- 1 ( t )  [A (t ) A ( t )  -  A (t) ] n  (t)  +  A"1 (t)B(t)u(t)  (54) 

From  (46)  of  the  previous  chapter, 

A  (t)  =  A (t)  A  (t )  -  A  (t )  J  (55) 

and  substitution  into  (54)  above  produces 

n (t)  -  Jn (t)  +  A-1 (t)B(t)u(t)  (56) 

where,  once  again,  J  is  the  constant  diagonal  matrix  of 
Poincare  exponents.  Equation  (56)  can  be  expanded  to 


hi  (t) 

Aihi (t) 

n2  (t) 

= 

A  j  h  2  ( t) 

+  A“l  (t )  B  (t)  u  ( t ) 

(57) 

• 

• 

n  (t) 
Ln 

• 

• 

A  n  (t) 
L  n  n  J 
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Since  the  system  allowed  diagonalization ,  which  in  this 
case  can  also  be  viewed  as  modal  separation,  it  became  ob¬ 
vious  that  feeding  back  a  portion  of  the  single  unstable 
mode  was  possible.  More  importantly,  the  diagonal vzation 
provided  for  a  relatively  easy  implementation  of  the  con¬ 
trol. 

Control  Law  and  Gain  Calculation 

As  stated  in  Chapter  III,  the  Poincare  exponent  anal¬ 
ysis  produced  two  stable  modes  and  one  unstable  planar 
mode.  Since  all  attracting  bodies  (earth,  moon,  sun)  were 
initially  assumed  to  lie  in  a  plane,  no  forces  which  would 
cause  out-of -plane  satellite  motions  were  present.  There¬ 
fore,  the  derivation  of  the  controller  was  based  on  the 

planar  modes  only,  i.e.  the  state  vector  for  controller 

T 

design  was  [R  p  R  p  ]  .  Since  the  controller  could  not 

x  x  y  y 

realistically  produce  instant anecus  position  changes,  the 

controller  was  physically  constrained  to  create  inertial 

velocity  changes  only.  The  difficulty  encountered  by  Smith 

m  his  work,  coupled  with  the  preceding  physical  con- 

strairt,  led  to  the  choice  of  the  simplest  possible  B  ma- 

T 

trix.  The  chosen  B  matrix  was  [0101]  .  Again,  since  in 

the  chosen  coordinate  frame  the  momenta  are  inertial  velo¬ 
cities,  this  choice  for  B  was  possible.  The  modal  equations, 
then,  were 
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Hi  (t) 

^  1 

Hi  (t) 

0 

• 

[0] 

n2  (t) 

A  2 

n2  (t) 

l 

= 

+A~ 1 

n3  (t) 

^  3 

n3  (t) 

0 

(t]_ 

>-< 

o 

1 

_n  4  (t) 

l 

u  (58) 


where  u  is  the  scalar  control,  also  chosen  for  simplicity. 
Choosing  a  scalar  control  and  a  B  matrix  as  defined  above 
is  equivalent  to  assuming  t'.u.t  the  control  is  applied  to 
the  satellite  at  45°  angles  from  the  rotating  +x  and  +y 
axes  shown  in  Figure  3  (since  p^  and  p^  receive  equal 
amounts  of  control) .  Defining  n3  as  the  unstable  component 
of  the  modal  vector,  the  control  law  became 


u  =  kn  3  (59) 

where  k  is  the  constant  control  gain.  Recalling  the  defi¬ 
nition  of  the  modal  variables,  (50)  ,  the  control  law  in 
physical  variables  was 

u  =  kAl1 (t) Sx ( t )  (60) 

where  A3x(t)  is  defined  as  the  third  row  of  A_1 (t) . 

To  calculate  a  value  for  k  which  will  theoretically 
stabilize  the  linearized  system,  (49)  ,  requires  a  solution 
of  the  unstable  component  of  (58). 

n3(t)  -  l3n3(t)  +  [A7J(t)  +  Al i ( t) ] kn 3 (t)  (61) 

(n3(t)  and  n*  (t)  were  defined  as  the  unstable,  planar  mode, 
while  m (t)  and  n2  (t)  constituted  the  stable,  planar  mode.) 
The  terms  A32(t)  and  A3i(t)  were  a  result  of  multiplying 
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AJ1(t)  by  the  D  vector  (the  32  subscript  denotes  third  row, 
second  column).  Rearranging  (61)  yields 


n3(t)  =  { A 3  +  [A 3 1  (t)  +  a3J (t) ] k } n 3  (t) 


(62) 


a  linear,  first  order,  time-varying  differential  equation 
for  n 3 (t) .  The  above  equation  can  be  solved  via  the  method 
of  an  integrating  factor.  For  notational  simplicity,  de¬ 
fine 

Y (t)  =  Ali  (t)  +  Ali  (t)  (63) 


By  the  method  of  an  integrating  factor 


From  the  Fourier  analysis  of  A_1  (described  in  the  next  sec¬ 
tion)  ,  it  was  discovered  that  A32<t)  and  A3£(t)  had  both 
constant  and  periodic  components;  therefore,  y(t)  can  be 
separated  into  a  periodic  function  and  a  constant. 


Y  (t )  =  yo  +  Y  (t)  (65) 

where  y  is  the  constant  and  Y(t)  is  the  periodic  function, 

o  p 

The  integrating  factor  expands  to 


e— / l A 3  +  kYQ  +  kYp  (t) ] dt 


(66) 


which  simplifies  to 

[e'{Xi  +  kYo)t] [e_/kYp(t)dt) 


(67) 
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Making  the  above  substitution  and  realizing  that  the  multi¬ 
plication  by  the  integrating  factor  produces  a  perfect  dif¬ 
ferential,  one  can  write  the  differential 


+  kVtn«--nVt)dtn  =  0 

and  the  solution  to  the  above  equation  is 

n,(t)  =  c[e'kVt,dtne<X’  +  kVh 


(68) 


(69) 


where  c  is  the  constant  of  integration.  The  stability  in¬ 
formation,  then,  is  contained  in  the  last  exponential  of 
(69) ,  since  the  periodic  term  will  merely  contribute  oscil¬ 
latory  motion  to  the  system.  The  shift  in  the  root  of  n 3 , 
then,  will  be  governed  by 


A3  =  A3  +  kyQ  (70) 

where  AJ  is  defined  as  the  new  root  location  after  control 
is  applied. 


_  Implementation  of  the  Controller 

The  final  implementation  of  the  controller  was  accom¬ 
plished  by  integrating  the  full  non-linear  equations  of  mo¬ 
tion,  described  by  (12)  and  augmenting  the  equations  with 
the  controller  described  in  (60).  Since  the  controller  was 
designed  for  the  linearized  system,  no  assurance  of  a  sta¬ 
bilizing  system  existed  when  the  augmented  non-linear  equa¬ 
tions  were  integrated.  The  vector  6x(t)  was  formed  by 
evaluating  the  state  of  the  non-linear  equations  at  time, 
t,  and  subtracting  the  state  of  the  periodic  reference 
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orbit  at  the  same  time,  e.g. 

6x(t)  =  x  ( t )  -  xref<t)  (71) 

The  matrix,  A  1 (t) ,  required  by  the  control  law,  (60), 
was  produced  as  foliows.  The  eigenvectors  of  the  system 
were  computed  by  integrating  (48)  with  a  small,  constant 
integration  step.  A  Fourier  series  was  then  constructed  to 
represent  the  elements  of  the  eigenvectors  at  each  inte¬ 
gration  step  such  that,  given  the  time,  one  could  produce 
the  eigenvectors  for  the  system  at  any  point  along  the 
orbit  from  the  Fourier  series  (Ref  2:60).  Since  A(t)  was 
composed  of  the  system  eigenvectors  in  accordance  with  (44), 
the  A (t)  matrix  was  then  produced  from  the  Fourier  series 
and  its  inverse  computed  at  many  points  evenly  spaced  along 
the  orbit.  The  components  of  the  A  1 (t)  matrix  were  then 
Fourier  analyzed  to  allow  the  A  1 (t)  matrix  to  be  construc¬ 
ted  at  any  given  time. 
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V .  Results,  Conclusions,  and  Recommendations 

Calculation  of  Controller  Gain 

Smith,  in  his  work,  derived  a  method  for  plotting 
what  he  termed  a  "quasi  root  locus"  for  a  variational  sys¬ 
tem  forequations  with  periodic  coefficients  (Ref.  8:29). 

If  the  coefficient  matrix,  A(t)  (also  called  the  plant  ma¬ 
trix  in  control  literature) ,  was  a  constant  matrix,  conven¬ 
tional  root  locus  techniques  could  be  employed  to  determine 
the  shift  in  system  roots  as  control  gain  is  varied.  In 
the  system  described  by  (49),  however,  the  matrix  A(t),  is 
time  periodic;  therefore,  a  conventional  root  locus  plot  of 
gain  versus  frequency  is  not  possible.  Equation  (49) 
though,  can  also  be  solved  by  Floquet  theory.  Smith  rea¬ 
lized  that  since  Poincare  exponents  are  determinants  of 
system  stability  in  a  Floquet  analysis,  solving  the  augmen¬ 
ted  state  equations  for  the  Poincare  exponents  for  various 
values  of  gain  results  in  a  "quasi  root  locus".  Implemen¬ 
ting  Smith's  theory  enabled  :he  calculations  required  to 
plot  the  root  locus  for  the  system  with  modal  control,  and 
the  plot  is  shown  in  Figure  5.  As  predicted  by  the  modal 
control  theory  of  Chapter  IV,  the  unstable  root  migrated 
to  the  left  on  the  real  axis  and  became  more  negative  as 
gain  was  increased.  Since  none  of  the  other  system  roots 
were  affected,  it  is  obvious  why  modal  control  can  be 
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Figure  5.  Quasi  Root  Locus 


viewed  as  a  pole  placement  technique.  The  values  of  yain 
used  were  k  =  0.0,  0.259,  0.359,  0.459,  0.559,  and  1.0.  As 
calculated  from  (70)  ,  the  minimum  gain  to  stabilise  the 
system  was,  in  fact,  k  =  0.259.  The  root  locus  analysis 
also  enables  one  to  plot  the  unstable  exponent,  A3,  versus 
control  gain,  and  this  plot  is  shown  in  Figure  6.  The 
linear  relationship  between  frequency  and  gain  was  unpre¬ 
dicted,  but  was  later  explained  by  the  gain  calculation  of 
(70).  Since  the  constant  portions  of  A32<t)  and  A-,il(t) 
were  found  to  be  -0.9857  and  -8.256,  respectively,  the  pre¬ 
dicted  slope  of  the  A3  vs.  gain  plot  of  Figure  6  was  9.2417. 
The  actual  slope  of  Figure  6  was,  in  fact,  9.2414,  thus 
verifying  the  validity  of  Equation  (70) . 

Evaluation  of  the  Controller 

As  stated  in  Chapter  IV  of  this  report,  the  control¬ 
ler  design  was  based  on  a  linearized  set  of  equations  and 
yet  the  controller  was  incorporated  into  the  non-linear 
equations  of  motion.  This  evaluation  technique  was  imple¬ 
mented  to  test  the  controller  in  the  most  realistic  model 
(available  to  this  study)  of  the  actual  dynamics  of  the 
system.  Care  was  exercised  in  the  choice  of  initial  condi¬ 
tions  to  insure  that  the  system  was  initially  inside  the 
region  where  the  linearization  was  valid. 

Recalling  the  form  of  the  modal  variables, 

n  (t)  =  A- 1  (t )  6x  (t )  (72) 


37 


/  .2  .3  .4  5  .6  7  Q  9  1.0 

GA  IN 


Figure  6.  X 3  vs  Gain 
3C 


The  ideal  controller  should  drive  the  modal  variables  to 
zero,  since  one  desires 


6x(t)  =  x ( t )  -  *ref(t)  =  0  (73) 

However,  since  only  ri  3  was  controlled,  one  would  anticipate 
the  stable  planar  mode,  m,  and  r)2»  would  oscillate  (due  to 
purely  complex  roots)  while  the  negative  real  portion  of 
the  unstable  mode,  rm,  exponentially  decreases.  In  the  ab¬ 
sence  of  control,  n.t  should  exponentially  increase.  These 
concepts  are  mathematically  demonstrated  by  recalling  (58) 
from  Chapter  IV  of  this  report.  If  u  =  0  (no  control) , 

(58)  reduces  to 


n  (t) 


A  i 


A  2 


[0] 


[0] 


A  3 


n  (t) 


(74) 


which  is  a  constant  coefficient  set  of  equations.  Since  A3 
and  A 2  are  purely  complex  conjugates  and  A3  and  A.,  are  a 
real  conjugate  pair,  the  theoretical  solution  to  (74)  is  of 
the  form 


ni  (tT 

HioCostAit)  +  ri2osin(A1t) 

12  (t) 

n2ocos(A1t)  *  niosinUjt) 

n a  (t) 

A  3  t 

n  3  o© 

_h4  (t^ 

—  A  3 1 

_ruoe  3  _ 
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where  i)io;  'i2o  ,  030,  and  n /.  u  are  the  initial  conditions  de¬ 
fined  by  (72)  at  t  =  0.  Figures  7  and  8  are  plots  of  the¬ 
oretical  and  computed  magnitudes  of  the  modal  variables  for 
two  orbital  periods  when  no  control  is  present  (i.e.  k  = 
0.0).  The  initial  conditions  (deviations  from  the  refer¬ 
ence  periodic  orbit),  6x  ( 0) ,  were  arbitrarily  chosen  to  be 
-1.496  km  in  x  and  y  position  and  -0.00298  cm/sec  in  x  and 
y  momenta.  Recalling  the  orbital  period  was  0.508  in  the 
units  of  Chapter  [I.  Figure  7  and  8  plot  two  orbits  (two 
months).  As  shown  by  the  two  figures,  the  computed  solu¬ 
tions  track  the  theoretical  solutions  quite  well. 

When  the  control  was  applied,  the  same  behavior  was 
seen  in  the  uncontrolled  moce,  aj  expected,  but  the  con¬ 
trolled  mode  was  forced  to  an  approximately  exponential 
decay.  Figures  9  and  10  are  plots  of  the  modes  with  the 
controller  activated.  The  plots  cover  four  orbits  about  L3 
and  the  controller  gain  was  k  =  0.359.  Initial  conditions 
chosen  for  this  case  were  -149.6  km  in  x  and  y  position, 

.  and  -2.978  cm/sec  in  x  and  y  momenta.  To  demonstrate  the 
effect  of  increasing  the  gain  of  the  controller.  Figures 
11  and  12  show  plots  with  identical  initial  displacements 
as  those  in  Figures  9  and  10,  but  the  controller  gain  was 
increased  to  K  =  0.559. 

The  oscillations  around  a  pure  cosine  or  sine  curve 
of  hi,  and  n 2 ,  and  the  overshoot  past  the  zero  level  of  ru , 
may  both,  in  part,  be  attributed  to  the  mathematical  con¬ 
struction  of  the  controller.  The  applied  modal  control 
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Figure  8.  Theoretical  and  Computed  Unstable  Mode  (2  orbits) 
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docs,  in  fact,  change  only  the  A3  frequency  as  seen  by  the 
root  locus  of  Figure  5.  Subsequent  analysis  illustrates 
that  all  system  eigenvectors  change,  and  this  change  in  the 
eigenvectors  causes  a  deviation  in  the  modal  behavior. 

This  deviation,  however,  is  insignificant  as  long  as  the 
system  is  stabilized.  Figures  9-12  do  portray  a  stabilized 
system. 

When  longer  periods  of  controller  operation  were  at¬ 
tempted,  undamped  oscillations  began  to  occur  in  aJ i  modes. 
Figures  13  and  14  portray  such  a  case.  The  magnitudes  of 
the  initial  conditions  chosen  for  this  case  were  much  smal¬ 
ler  than  those  of  Figures  9-12  (-1.496  km  in  x  position, 
+1.496  km  in  y  position,  +0.00298  cm/sec  in  x  momentum  and 
-0.02978  cm/sec  in  y  momentum  for  this  case),  and  the  chosen 
gain  was  k  =  0.359.  Figures  13  and  14  represent  two  years, 
or  24  orbits,  of  operation.  In  Figure  14,  the  oscillations 
occur  in  both  n3  and  n<*  and  the  oscillations  appear  to  be 
sinusoidal  with  a  linearly  increasing  magnitude.  It  is 
-  also  noteworthy  that  the  controller  did  overcome  the  ini¬ 
tial  displacements,  but  then  started  to  oscillate  after  at¬ 
taining  the  desired  zero  level  for  n 3  and  ru, . 

To  isolate  the  cause  of  the  oscillations,  the  same 
case  was  rui  (k  =  0.359,  24  orbits),  but  the  initial  dis¬ 
placements  were  removed  such  that  6x(0)  =  0.  Figures  15 
and  16  demonstrate  that  the  oscillations  were  actually 
present,  though  small  (note  the  scale) ,  before  they  were 
apparent  in  Figures  13  and  14. 
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Several  possibilities  were  postulated  for  the  oscil¬ 
lations  seen  in  the  results  thus  far.  The  first  of  these 


was  an  "out  of  sync"  condition  between  the  reference  orbit 
of  the  satellite  and  the  integrated  non-linear  equations  of 
motion  augmented  with  the  controller.  Since  the  L3  orbit 
is  represented  by  a  Fourier  series,  it  is  possible  that  the 
orbital  frequency  and  the  frequency  of  the  integrated  equa¬ 
tions  might  become  mismatched.  To  illustrate,  assume  that 
the  integrated  frequency,  u>,  is  mismatched  by  a  small  fac¬ 
tor,  e,  such  that  the  periodic  motion  is  governed  by  a  term 
such  as 

sin(ot  +  et)  (76) 

Expanding  (76)  yields 

sin  (oat) cos  (et)  +  cos  (wt)  sin  (et)  (77) 

Since  e  is  small,  cos(et)  s  l  and  sin(et)  s  et.  Substitu¬ 
ting  these  approximations  into  (77)  yields 

sin(mt)  +  cos(u>t)et  (78) 

which  would  produce  the  observed  periodic  motion  plus  an 
increasing  amplitude  term  (mixed  secular  term) .  To  inves¬ 
tigate  the  possibility  that  the  oscillations  were  a  result 
of  this  frequency  mismatch,  the  system  was  allowed  to  ope¬ 
rate  for  one  year  (12  orbits)  with  no  control.  Figures  17 
and  18  show  this  case  which  further  demonstrates  the  be¬ 
havior  shown  in  Figures  7  and  8.  Unlike  Figures  7  and  8, 


however,  Figures  17  and  18  represent  a  case  with  zero  ini¬ 
tial  conditions.  The  modes,  then,  were  well  behaved  for 
approximately  seven  months  at  which  time  the  unstable  eigen¬ 
value,  eA:jt,  had  increased  enough  to  dominate  system  be¬ 
havior.  Figures  19  and  20  further  illustrate  the  need  for  ’■ 

j 

control.  After  two  years  with  no  control,  even  though  the  ] 

1 

initial  conditions  were  zero,  the  system  is  totally  out  of 
the  linear  regime  and  the  calculations  were  nonsensical. 

Another  possibility  for  the  oscillations  that  merited 
investigation  was  the  size  of  the  integration  step.  If  the 
step  size  were  too  large,  the  integration  might  be  compu¬ 
ting  an  inaccurate  state  vector,  which  then  is  used  to  cal¬ 
culate  6x(t).  When  the  step  size  was  halved  for  the  case 
shown  in  Figures  21  and  22,  the  oscillations  remained. 

This  case  also  had  the  controller  inoperative  and  6x(0)  =  Cf. 

Inclusion  of  Eccentricity  and  Inclination  Effects 

Chapter  II  of  this  report  described  the  periodic 
orbit  for  the  moon  which  was  constructed  by  Wiesel.  Re¬ 
sults  presented  thus  far  are  based  on  this  lunar  model. 

The  Wiesel  lunar  Oibit  program,  though,  is  also  capable  of 
providing  the  eccentricity  and  inclination  effects  of  the 
moon.  With  these  two  forcing  functions  included  in  the 
lunar  model,  the  controller  was  truly  tested  on  a  very  re¬ 
alistic  model  of  the  earth-moon-sun  system.  With  a  control 
gain  of  0.459  and  the  initial  conditions  set  to  zero,  the 
system  was  operated  for  four  months  (4  orbits) .  Figures  23 
and  24  depict  the  behavior  of  the  modes  for  this  case. 
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jure  21.  Stable  Mode  (4  orbits,  6x(0)  =  0,  smaller 
gration  step  size) 
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When  Figures  23  and  24  are  compared  with  Figures  15  and  16, 
one  notices  the  oscillations  are  much  larger  in  Figures  23 
and  24,  and  the  oscillations  of  n3  and  n4  appear  to  have  a 
frequency  approximately  equal  to  the  orbital  period. 

The  four  components  of  <5x(t)  which  are  used  to  con¬ 
struct  n  are  also  plotted.  Figure  25  portrays  Sx*  (t)  and 
<$x2  (t)  while  Figure  26  shows  6x3  (t)  and  6x4  (t)  .  The  curves 
appear  to  reflect  the  same  type  of  behavior  demonstrated  by 
the  system  prior  to  the  inclusion  of  the  forcing  terms. 

The  controller  does  maintain  the  6x(t)  components  around 
the  zero  level,  but  a  periodic  forcing  term  is  readily  ap¬ 
parent.  For  the  four  orbits,  oscillations  observed  were 
approximately  4400  km  in  x  position  and  6000  km  in  y  posi¬ 
tion. 

without  the  lunar  inclination  and  eccentricity  ef¬ 
fects  included,  a  calculated  control  cost  is  superfluous 
since  it  is  not  based  on  this  study's  most  realistic  model. 
Control  costs  were  computed  for  the  case  described  above 
which  is  shown  in  Figures  23-26.  The  control  cost  was  com¬ 
puted  by  calculating  the  magnitude  of  the  control  applied. 
Since  control  is  actually  an  acceleration  (A  velocity) , 
t 

Av  =/  |control  acceleration j dt  (79) 

0 


where  Av  is  the  change  in  velocity.  To  compute  an  average 
control  acceleration,  then, 


_  1  ^ 

a  *  —  /  | control  acceleration | dt 

*  0 

where  a  is  the  average  acceleration. 


(80) 
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Figure  26.  6xj  and  fix*.  Lunar  Inclination  and  Eccentricity 
Effects  Included  (4  orbits) 
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This  equation  was  then  integrated  along  with  the  other  six 
equations  of  motion  to  form  an  integrated  control  cost.  At 
the  end  of  the  four  month  period  of  operation,  the  control 
cost  for  the  case  described  by  Figures  24-27  was  1.026  X 
10_6  m/sec2  which  is,  after  conversion,  an  average  accelera¬ 
tion  of  1.047  X  10  7  g.  By  comparison,  Farquher  (Ref  3:63) 
cites  earth-synchronous  satellite  control  costs  of  1.51  X 
10“7  g. 

Conclusions 

As  evidenced  by  Figures  9  and  10,  the  controller  does 
exhibit  predicted  theoretical  behavior  by  driving  n3  and  i|4 
to  zero  and  not  affecting  m  or  112-  The  oscillations  pre¬ 
sent  in  the  behavior  of  the  modes  for  longer  runs  are  not, 
therefore,  induced  by  the  controller.  Two  possibilities 
have  been  postulated  for  the  oscillatory  modal  behavior. 

As  previously  stated,  since  the  reference  periodic  orbit 
about  L3  was  constructed  in  single  precision  arithmetic 
with  an  integrator  accurate  to  10  1 0 ,  and  then  represented 
by  a  Fourier  series  (an  approximation  technique  when  the 
Fourier  series  is  truncated) ,  it  is  possible  that  the  ref¬ 
erence  orbit  inaccuracies  caused  the  reference  orbit  to  be 
out-of-sync  with  the  actual  integrated  and  controlled 
orbit.  The  other  possibility  is  that  the  moon's  orbit  cre¬ 
ates  the  oscillations  seen  on  the  plots.  In  all  cases,  the 
period  of  the  oscillations  was  approximately  the  L3  orbital 
period  which  is  also  the  period  of  the  moon's  orbit.  In 
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the  case  of  the  plots  of  and  02  for  the  different  cases, 
it  appears  that  the  periodic  forcing  function  of  the  moon  is 
modulating  the  underlying  predicted  (and  observed)  sine 
wave.  For  1)3  aid  0  ,, ,  the  same  forcing  function  appears 
superimposed  on  the  predicted  exponentials.  if  the  lunar 
orbit  does,  in  fact,  cause  the  oscillations ,  the  oscilla¬ 
tions  can  be  viewed  as  a  result  of  resonance  between  the 
lunar  and  L3  orbits. 

The  control  cost  cited  for  the  case  shown  by  Figures 
23-26  was  31  percent  lower  than  the  costs  incurred  by 
satellites  in  synchronous  orbits.  Since  control  costs  were 
the  primary  yardstick  for  feasible  of  modal  control,  one 
concludes  that  modal  control  is  a  viable  scheme.  The  orbit 
about  L3  provides  an  excellent  location  for  sensor-bearing 
satellites,  particularly  when  coverage  of  the  polar  regions 
is  of  paramount  importance.  One  can  also  state  that  the 
same  control  technique  could  be  applied  to  satellites  in 
orbits  about  L4  and  L5,  providing  the  desired  total  global 
coverage . 

Recommendations 

Follow-on  studies  should  construct  the  periodic  ref¬ 
erence  orbit  about  b3  with  double  precision  arithmetic  to 
insure  the  out-of-sync  condition  mentioned  previously  does 
not  exist  or  is  minimized.  The  resonance  condition  (if  it 
exists)  is  a  function  of  the  dynamics  of  the  problem  and 
therefore  must  be  tolerated. 
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A  further  attempt  to  alleviate  the  oscillation  (ring 
ing)  problem  would  be  to  use  some  sort  of  state  estimation 
technique  rather  than  the  purely  deterministic  dynamics  of 
this  study.  State  estimation  should  cure  any  frequency 
mismatch  problems  that  might  exist  and  provide  an  easy  ve¬ 
hicle  for  optimizing  the  control  law.  State  estimation 
techniques  have  been  used  extensively  in  the  studies  of 
control  of  large  space  structures  and  the  theory  should 
directly  apply  to  this  problem. 

Since  the  Poincare  exponents  of  one  planar  mode  and 
the  out-of-piane  mode  are  purely  imaginary,  one  can  only 
say  that  they  are  marginally  stable.  Control  theory  pre¬ 
dicts  pure  oscillation  for  roots  on  the  imaginary  axis  and 
the  expected  oscillatory  motion  for  these  modes  does,  in 
fact,  occur.  Future  studies  should  investigate  adding 
damping  to  these  modes  to  move  the  roots  off  the  imaginary 
axis  and  into  the  left-half  of  the  complex  plane,  thereby 
producing  stable  modes  rather  than  marginally  stable  modes 
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